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Abstract 



We consider a finite number of coupled oscillators as an adaptation of the Kuramoto model of 
populations of oscillators. The synchronized solutions are characterized by an integer m, the 
winding number, and a second integer /. Synchronized solutions of type (m, / = 0) are all stable, 
and an explicit perturbative expression for these for large values of the coupling constant K is 
presented. For low K, these solutions appear at certain specific values, each merging with a solution 
of type (m, / = 1), both these solutions being stable for K close to the relevant value. The (m, 0) 
solution continues to be stable for larger K, while the (m, 1) solution, on continuation in K becomes 
unstable and then merges with a new type (with a different m and /) of unstable solution. The (m,0) 
type solutions for large K are in the nature of phase waves traveling round the ring. All the stable 
synchronized solutions are entrained by an external periodic driving, provided that the driving 
frequency is sufficiently close to the frequency of the synchronized population. A perturbative 
approach is outlined for the construction of the entrained solutions. For a given amplitude of 
driving, there is a certain maximum detuning between the two up to which the entrained solution 
persists. The question of stability of the entrained solution is addressed. The simplest situation 
involving three oscillators is investigated in details where the onset of synchronization is seen to 
occur through a tangent bifurcation at some critical value of K. Immediately before the appearance 
of the first synchronized solution, the system exhibits intermittent chaos with a typical scaling for 
the duration of the laminar phase. Under a periodic driving with an appropriately limited detuning, 
there occurs entrainment of the chaotic solution. 

PACS numbers: 05.45.Xt, 87.10. +e 

Keywords: nonlinear oscillator, synchronization, circadian rhythm, entrainment, Kuramoto model 



*Electronic address: tarun.roy@saha.ac.in 
^Electronic address: la'l@vsnl.coml 



2 



I. INTRODUCTION 



Synchronization in coupled non-linear oscillators has been observed and investigated for 
a long time now. Starting from the observation by Christian Huyghens in 1665 of 
synchronization in weakly coupled pendulum clocks, a wide variety of phenomena involving 
synchronization has been explored, including the spectacular synchronization of the flashing 



of fireflies jsl- In recent decades, a theoretical understanding of synchronization has played a 
definitive role in discerning and explaining the mechanism underlying biological rhythms of 
a wide variety. Thus, populations of self-oscillatory neuronal cells with mutual interactions 
exhibit the phenomenon of phase locking, where all the cells oscillate in phase, as in the 
cardiac pacemaker, studies on which were initiated in a pioneering work by Peskin 4|. 
Similar synchronization phenomena are found in the insulin-secreting cells in the pancreas, 
in populations of cells controlling such activities as breathing, walking, and running, and in 



the rhythmic activites of the intestines 



A simple but highly fruitful model of an interacting population of nonlinear oscillators, pro- 
posed by Kuramoto (see, e.g., |3]), opened up a new era in the analysis of large populations 
of interacting nonlinear oscillators where synchronization appears as a macroscopic mode of 
behaviour for such a population. The onset of synchronization could be studied in the model 
in terms of a macroscopic order parameter relating to the population. A vast literature is 
now availa ; 
Winfree js, 

diverse areas in the physical and biological sciences on the basis of the model. A comprehen- 
sive review of various aspects of the Kuramoto model with applications in neural networks, 
Josephson junction arrays, laser arrays, charge density waves, and laser oscillators, is to be 



e on the Kuramoto model (which followed a class of models proposed earlier by 
qJ) and its variants, and on the explanation of synchronization phenomena in 



found in 



iq|. 



Another area of study where the Kuramoto model is potentially of great use is the phe- 
nomenon of entrained synchronization where the effect of externally imposed periodic driv- 
ing on a synchronized population of coupled nonlinear oscillators is explored with the help 
of the model. Entrainment is a term generally used to denote the influence of two or more 
oscillating systems on one another whereby they oscillate at a common frequency. According 
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to this definition, the phenomenon of synchronization itself is an instance of mutual entrain- 
ment of a population of oscillators. However, we shall use the term in a more specialized 
sense here, where it is meant to describe the modofication of the frequency of synchronized 
oscillations of a population of oscillators by an externally imposed periodic influence such 
that the latter makes the entire population oscillate with its own frequency. 

This phenomenon of entrainment operates in the circadian time-keeping system of organisms 
where the external influence (commonly refrred to as a zeitgeber in the biological literature) 
is the day-night cycle of light intensity that affects the endogenous oscillations by means of 
signals generated through appropriate optical receptors. 

The problem of exogenous entrainment of an endogenously oscillating synchronized popula- 
tions of oscillators described by the Kuramoto model has been considered by Sakaguchi 
and more recently by Ott et al [12]. The organization of the circadian timing system (CTS) 



and a simplified model of this syetem, has been described in [13|, where the pacemaker of 
the suprachiasmatic nucleus (SCN), a central component of the CTS, has been treated as a 
large but finite population of cells described by the Kuramoto model. 

While the approach initiated by Kuramoto and followed in a large number of subsequent 
investigations considers the thermodynamic limit N ^ oo, where stands for the number 
of interacting oscillators, and then looks for synchronized state in terms of a macroscopic 
order parameter, the behaviour of the population for a finite N has also been investigated, 
with results analogous to the thermodynamic approach, though exact results for finite 



are few in number. Maistrenko et al 1^ have considered the scenario just before the onset 
of synchronization with global coupling (i.e., each oscillator coupled to all the remaining 
oscillators in the population) for finite A^. The case = 3 is a well-known torus flow 
while N = 5 gives important insights into the phase flow leading to synchronization. In 



a subsequent work, Popovych et al [15[ demonstrated the existence of a high dimensional 
chaotic regime in the finite- A^ Kuramoto model, such high dimensional chaos being a common 
feature in coupled oscillator systems. 

Synchronization in finite linear chains of oscillators has been investigated in a number of 



papers (see 



16| and references therein) with nearest neighbour coupling, multiple neighbour 
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coupling, and also with global coupling. Travelling wave solutions were obtained in these 
models, where the range of coupling was found to affect the phase lags between the successive 
oscillators. Such chains of oscillators are relevant as models of central pattern generators 



(CPG) in certain organisms (see 17| and references therein). 



In the present paper we consider a finite number of coupled oscillators with a ring-like 
connectivity, based on a nearest-neighbour interaction, where each oscillator is represented 
by a time-dependent phase as in the Kuramoto model. Rings of coupled oscillators were 
considered in {g], where each oscillator was an weakly attracting limit cycle (other possible 
degrees of freedom having been eliminated by a centre-manifold reduction). Phase wave 
solutions travelling round the ring were found, and a stabilty analysis of these phase waves 
was also performed. Oscillating systems coupled into a ring-like topology serve as a model 
for a wide array of applications such as locimotion gaits, rings of semiconductor lasers, and 
circular antenna arrays [l^. Rogge and Aeyels 19 1 considered a chain of oscillators with 
unidirectional nearest neighbour couplings and derived an algorithm for the construction all 
the phase-locked solutions of the system. They also made a detailed stability analysis of 
these solutions, deriving interesting results. 



Our model is similar to that in 



19| with, however, a bidirectional coupling as in [6] , where 



the system consists of oscillators as in the Kuramoto model, each represented by a time 



201 ]. who derived 



dependent phase. The same model has been considered by Ochab et al 
an algorithm for the synchronized solutions of the system. We follow the approach of these 
authors and work out the stable solutions in the limit of large K, where K stands for 
the coupling constant. For sufficiently small values of K, the system does not admit of 
synchronized solutions. At larger values of K, there appear synchronized solutions where 
each such solution is characterized by two integers m, /(!</< N). The number and 
types of synchronized solutions change as K is made to assume increasing values. This is 
illustrated in details for the simplest non-trivial case A = 3, where a stable and unstable 
solution are seen to appear through a tangent bifurcation at a certain critical value of K 
(= Kc), the two solutions being characterized by the integers m = 0, / = 1. As if is made to 
increase, the stable (0, 1) solution changes to a (0, 0) type solution at a certain vealue of K, 
while the unstable (0, 1) solution changes into a (m = 1, / = 2) type solution. Similar other 
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scenarios, for other (m, I) type solutions, appear at still larger values of K. For K slightly 
below Kc, the system exhibits intermittent chaos, characterized by alternating laminar and 
turbulent phases, with a characteristic scaling of the duration time of the laminar phase and 
the Lyapunov exponent. 

An explicit perturbative construction of solutions is possible for sufficiently large values of 
K, where it is seen that only the / = type solutions (with various different values of m; 
for = 3, the only possible value of m for / = is m = 0) are stable. The / = solutions 
for various values of m appear as stable phase waves traveling round the ring of oscillators. 
While our results constitute an extension of those in 20|], a complete picture describing all 



the stable and unstable synchronized solutions for all values of K, and for arbitrary N, is 
yet to emerge. 

Finally, we extend the model to include an external periodic driving on the system and 
study the entrainment properties as a function of the apmlitude (A) and frequency (Q) 
of driving. For any given frequency Q, sufficiently close to the free-running synchronized 
frequency (u) of the system, there exists a critical amplitude Aq above which entrainment 
to the external frequency Q occurs. We present a perturbative construction of the entrained 
solutions starting from the synchronized solutions obtained in our model. It is also found 
that the question of stability of the entrained solutions reduces to simple terms since stability 
is seen to be determined by a constant rather than a periodic matrix. 

The plan of the paper is as follows. In sections [Til and IIIII we formulate the Kuramoto 
ring problem consisting of N oscillators with bidirectional nearest neighbour interaction, 
and set up the basic equation describing the synchronized solutions. Section [IV] includes the 
formulation of the stability problem for the possible synchronized solutions of the system. An 
explicit construction of solutions of type (m, / = 0) is presented in section |V] for large values 
of the coupling constant K, where these solutions appear as stable phase waves travelling 
round the ring. Section I VI I explores the low K regime where the synchronized solutions 
with various values of m and / appear in succession. Detailed considerations for the case 
A^ = 3 are presented in section lylTl for a particular choice of the frequencies Ui [i = 1, 2, 3) of 
the oscillators. Section IVIIII describes the phenomenon of intermittent chaos characterizing 
the system just before the appearance of the first synchronized solution as K is made to 
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increase from a low value. Finally, section |IX] addresses the problem of entrainment of the 
synchronized solutions under a periodic external driving including that of stability of the 
entrained solutions. The entrainment scenarios for small and large values of the driving 
amplitude A are briefly compared. The phenomenon of entrainment of intermittent chaos 
is also indicated on the basis of numerical results. 



II. RING OF OSCILLATORS WITH NEAREST NEIGHBOUR COUPLING 



We start with number of oscillators with time-dependent phases 6i {i = 1,2, . . . , N), 
described by equations 

K 

9i=uji + —{sin{e,+i-ei) + sin{9i_i-9i)), {i = 1,2, . . . , N), (la) 

where we have assumed nearest neighbour interactions and where the ring boundary condi- 
tion is expressed in the form 

(^N+i = di, 9q = 9n. (lb) 

In these equations, Ui are the natural frequencies of the oscillators making up the ring, 
and X is a coupling constant {K > 0) which will be treated as a parameter characterizing 
the model. Following (20], we derive in the next section the algorithm for obtaining the 
synchronized solutions of the system and then work out the stable solutions explicitly in the 
large K limit. 

We look for synchronized solutions in the model, of the form 

9^{t) = tut + (2) 

where (pi are time-independent phases, and u stands for the common frequency of the oscil- 
lators in the synchronized state. Substituting in eq. (ITal) . one obtains 



K 

u = Ui + —{sin{(j)i+i - (pi) + sin{(l)i_i - (f)i)), {i = 1,2, . . . , N). (3) 

Adding up all the equations in eq. ([2]), one finds that oj, the frequency of synchronized 
oscillation, has to be the mean of the oscillator frequencies: 

N 



^E'^- (4) 



N 

i=l 
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and the time-independent phases {i = 1,2, N) satisfy 

2 

{sin{(l)i+i - (j)i) + sin{(j)i^i - 0^)) = -^^i {i = 1,2, . . . , N). (5a) 

In this equation, Aj denotes the frequency deviation of the zth oscillator from the synchro- 
nized mean frequency, 

A,=u,-u iz = l,2,...,N), (5b) 
and the time-independent phases (pi satisfy the ring boundary conditions 

(pN+l = 01, 00 = 0Af- (5c) 

III. SYNCHRONIZED SOLUTIONS IN THE RING MODEL 

The basic equation (l5a|l determining a synchronized solution involves the time-independent 
phases (pi only through 

Pi = sin ipi, 'ipi = <pi- (z = 1,2, . . . , A^). (6a) 

in the form 

Pi+i -Pr = -^Ai, {t = l,2,...,N), (6b) 

where 

Pn+i=Pi- (6c) 



Evidently, the equations fl6b[) are not independent of one another since they add up 
to zero. Choosing p = pi (say) as a parameter for the time being, all the other p^'s are 
determined as 

2 

Pi=P+j^Si, {t = l,2,...,N), (7a) 

where 

i-l 

s, = -5^A,, {i = 2,...,N), si = 0, (7b) 
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and where, in equations (17a|) . the first equation results from the definitions of p and si. 



With all the piS so determined in terms of the parameters of the model (the coupling 
constant^ and the frequencies Ui) and of one can finally determine p by noting that the 
phase differences tfji have to satisfy 

N 

J]^i = 27rm, (8) 

where m is an integer, to be referred to as the winding number, on which more later. 
This equation involves pi through the relations pi = sin ipi (see equation f l6a|) ). and hence 
determine p by virtue of eq. ( !7a|) . 

We note that ipi is the inverse sine of Pi which may lie in either of the two intervals [— f , + f ) 

and [|, ^), i.e., with reference to the unit circle on which angles are depicted, either on the 

right half-circle or the left half-circle. For any given x in the range — l<z<+l,we denote 

these two values of the inverse sine of x as a{x) and I3{x) respectively. Here the end-point 

^ is not identified with — | since it may arise as a limiting value of an angle belonging to 

the left half circle. Thus, the equation (jH]) determining p (and hence a synchronized solution 

for the system) is finally seen to be of the form 

2 2 2 2 

{a{p + — SiJ + ■ ■ ■ + a(p + — Si^.J) + (/3(p + — s^J + ■■■ + [3{p+ ^^iJ) = 2m7r, (9a) 

where / is any integer ranging from to and where it is understood that, for / = 0, the 
second group of terms in the left hand side of the above equation is absent while, for / = A^, 
the first group of terms is absent. Here ii, . . . , in-i and ji, . . . , are two disjoint subsets 
of the set of integer suffixes 1, . . . , A^. A synchronized solution is thus characterized by the 
integers m and I, where / denotes the number of phase differences [ipi) that belong to the 
interval [f ,3|) (i.e., to the left half of the indicator circle). 

For the particular case of a (m,0) type solution, the above equation determining p assumes 
the form 

^ 2 

^a{p+ —Si) = 2mn. (9b) 

i=l 

We choose the number of oscillators A^ and the frequencies uji (and hence s,, i = 1, . . . ,N) 
as given parameters for a ring of oscillators, and treat the coupling constant as a variable 
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parameter, and we address the task of finding the synchronized solutions of the system 
for various possible values of K. Any such synchronized solution corresponds to specific 
values of the integers l{0 < I < N) and m and to a certain subset ji, . . .ji of the set 
1, . . . , (the complementary subset ii, ■ ■ ■ ,iN-i is automatically determined). For given m 
and / (1 < / < A^), there may be several different solutions of type (m, I) depending on 
which of the N phase differences ipi i = 1^ . . . , N are in the left half-circle and which in the 
right half-circle. 

For a given /, since the upper and lower limits of the phase differences are known, one can 
determine the range over which the integer m can vary, where one finds 

N I N I , ^ 

where the equality signs are not included since these can arise only for a degenerate solution 
with all the ^j's equal in the special case when all the oscillators have the same frequency. 

In particular, in the case of a (m, 0) type solution, m can be any interger in the range 

N N , , , 
<m< —. (10b) 

The integer m is referred to as the winding number characterizing the solution since, starting 
from the phase for any particular oscillator, say, 0jv and following the rotations on the unit 
circle by the successive phase differences i/^i, • • •, one comes back to the same oscillator 
after m windings of the circle. 



IV. STABILITY OF A SYNCHRONIZED SOLUTION 

For any given synchronized solution determined as above, a rigid rotation of all the phases 
0j (z = 1, . . . , A^) on the circle by an arbitrary angle $, corresponding to the transformation 
0j — )■ 0j + $, does not alter the phase differences ipi. In other words, a synchronized solution 
is actually a member of a one-parameter family of solutions defined by rigid rotations. This 
immediately goes to show that the synchronized solutions are all neutrally stable, i.e., the 
matrix (M) describing the hnearized equation obtained from eq. (fTa|) always possesses at 
least one zero eigenvalue. 
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Defining qi = cosipi, for the synchronized solution under consideration, the elements Mij (1 < 
hj < N) of this matrix are seen to be given by 



K K 

Ma = - — {qi + ^j+i), Mi^i+i = y^j+i, Mij 



{l<t,j<N), 



(11) 



where, once again, the indices are to be interpreted in accordance with the ring boundary 
conditions (A^ + 1 = 1, = A^), and where all elements not specified by (fTTl) are zero. For 
instance, with N = 5, the stability matrix is 



M 



K 
~2 



\ 



(?1 + ?2) 


(12 










-{<i2 + qs) 














94 










-(94 + gs) 










<l5 








(12) 



95 

-{q5 + qi)J 

It is easy to see that, for any A^, this matrix possesses a zero eigenvalue regradless of the 
values of gi, . . . , since, for the column x = (1,1,..., 1)""", one has Mx = 0. 



Moreover, in accordance with Greshgorin's theorem 19|, |2l|, the following results hold: (a) 
if all the ipi^s belong to the right half circle, then all eigenvalues of M are less than or 
equal to zero, in which case the synchronized solution under consideration is orbitally stable 
or orbitally neutrally stable, i.e., in other words the (m, 0) type solutions are, in general, 
orbitally stable; (b) if, on the other hand, all the i(ji^s are in the left half circle, the solution 
under consideration is orbitally unstable or, at most, orbitally neutrally stable, i.e., (m, 
/ = A^) type solutions are, in general, orbitally unstable. 

In addition, based on Greshgorin's theorem, the following result can be stated by employing 
a continuity argument: a (m, 1) type solution, merging with a (m,0) type solution at 
some particular value of K (see below), can remain orbitally stable as K is made to vary 
continuously up to a certain limiting value. 
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V. PHASE WAVES FOR LARGE K 



For sufficiently large K, for which all the terms -^Sj (z = 1, . . . , A^) are small and p is bounded 
away from ±1 for a solution of type (m, 0) (see below), one can make the expansion 

2 2 

sin~^{p+ —Si) ~ sin^^p-\ Sj, (13) 

K K^\~p^ 

where higher order terms in the expansion are ignored. 

Thus, from eq. (l9bl) the value of p for a (m,0) type solution at large K works out to 

. Imix 2 

p^5zn^-^^s„ (14) 

i=l 

where, once again, only terms of the order of are retained. This gives the result 

• -1/ 2 , 2m7r 2 , , ,2m7r, ,^ ^ , 

^i = sin {p+ —Si) = - —{s- Si)sec{^^), (15a) 

where 

AT 



i=l 

In the limit K — )■ oo, one therefore obtains the synchronized (m,0) type solution in the form 
of a one-parameter family 

2rmr , ^ 

9i = oot+ -j^i + 5, (16) 

where the parameter 5 stands for the constant part of the phase (independent of t and z), 
which can be chosen arbitrarily. Evidently, such a solution is in the nature of a phase wave 
running round the ring of oscillators, corresponding to a uniform phase advance between 
successive oscillators. Synchronized solutions for lower values of however, do not appear 
as phase waves since these are characterized by a non-uniform phase advance. 

Incidentally, the range of variation of m given by (llObI) shows that p (eq. (fT4l) ) is indeed 
bounded away from ±1 as we assumed to start with. 
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VI. SYNCHRONISED SOLUTIONS: THE LOW-i^ SCENARIO 



Continuing with the (m,0) type solutions for lower values of K, we will now see that a 
necessary condition for the existence of such a solution is that the coupling constant K is 
to be larger than a certain minimum value Kq. For K > Kq, a (m,0) type solution appears 
at a certain threshold value Kq"^\ which continues as a stable solution at larger values of K 
till, at a sufficiently large value of K the solution assumes the form of a phase wave given 
by (HSiD. 

Since the phase differences are given by sm-i/^j = p + ^Si, a necessary condition for the 
existence of a synchronized solution is that —1 < p + < 1 for alH = 1, . . . , A^. In other 
words, denoting by x and y the maximum and minimum values in the set {si, . . . , sat}, p 
and K must satisfy the conditions 

2 2 

P+ < 1, P + -^y > -1. (17) 



Thus, defining 



p+iK) = 1 - f a;, p.{K) = -1 - ^y, (18) 



the possible values of p and K for which synchronized solutions exist correspond to points 
lying in between the graphs for p^{K) against K, and p-{K) shown schematically in fig. [H 
where the point of intersection of the two graphs corresponds to 

X ~\~ y 

K = Ko = X - y, p = po = . (19) 

x-y 

As an example of how a (m, 0) type solution emerges as K is made to increase from the value 
Kq, we consider the emergence of the (m = 0,/ = 0) solution. This solution is determined 
by the relation 

N 

n(r) 4- 

K 



J]a(p+|;s,) = 0, (20) 



where, recall that a{x) stands for that value of sirT^x which lies in the range — f < < f 
(the inclusion of the end-point | does not involve any inconsistency since, as we see below, 
this value is relevant only at a point where a solution of a given type passes continuously 
into one of a different type). 
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FIG. 1: Graphs showing the variation of p^{K) and p-^-{K) with K (schematic); the two graphs 
intersect at p = po, K = Kq; the values of p and K for any synchronized solution has to correspond 
to a point lying in between these two curves; thus, for the existence of a synchronized solution, K 
has to be larger than Kq. 

Evaluating the left hand side of eq. (120|) with p = po and K = Kq , one obtains the expression 



If this expression is negative then, considering the two graphs of p^{K) and p^{K) referred 
to above, the equation ( l20l) can be satisfied only if (p, K) lies on the graph for p^(K) 
while, conversely, if P is positive then eq. (l20l) can be satisfied only on the graph of p-{K). 
This is because p^{K) increases monotonically with increasing K while p^{K) decreases 
monotonically. 

Let us assume, for the sake of illustration, that the frequencies of the oscillators are such that 
P < 0. Then, starting from the point (^0,^0), as we move along the upper graph, a point 
A(ir*^°\p_i_(ir*^°^)) will be encountered at which eq. ( |20l) will be satisfied. The value K = K^'^^ 
(with a slight change in notation, K^^^ is shown as Ki in fig. [1] for the sake of convenience) 
then corresponds to the emergence of the (0, 0) type solution since, for K < at least 




(21) 
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one of the expressions p + -^Sj {i = 1, . . . , N) will be greater than +1. On the other hand, 
for a K satisfying K > K^^\ one will have a (0,0) type solution with p{K) satisfying (|20|) . 
Here we have used the abbreviated notation K^^^ for where K = Kq"^^ corresponds to 
the appearance of a (m, 0) type solution. 

What is interesting to note in this connection is that, when K is slightly less than 

one can, by continuity, have a solution to eq. fl9a|) with 1 = 1, where now the largest of the 

phase differences -^j will be slightly larger than |. In other words, a (m = 0, / = 1) solution 

transforms to a (m = 0,/ = 0) solution at = K^^\ the latter being the value of K for 

which the largest of the phase differences for either of the two solutions assumes the value 
— 

2 • 

Of the two synchronized solutions that meet at K = K^^\ the (0,0) solution is orbitally 
stable by Greshgorin's theorem since all the g^'s in the linear stability matrix M (refer to 
sec. HVT) are positive. For K close to K'^'^\ the associated (0, 1) type solution is also stable 
by continuity, since the two sets of phase differences {tfji {i = 1, . . . , A^)} are close to each 
other. 

As K is made to increase away from K^^\ the corresponding values of p for the (0,0) 
type synchronized solution trace out a graph depicted schematically in fig. [1] (recall that 
this illustration assumes P < 0, where P is given by the expression fl2T|) ). which is drawn 
with a solid line to indicate that it is an orbitally stable solution. Starting from the point 
A (K = K^^\p = pj^^(K^^^)), the graph approaches asymptotically the value p = in 
accordance with iHM for X — )■ oo. 

As indicated above, the corresponding graph for the (0, 1) type solution is seen to start at 
the same point A. Considering for the sake of concreteness, the case of = 3 (see se dVIip . 
the graph depicting the variation of p with K for the (0, 1) solution is seen to include a 
point B where it shows a backbending. At the point B (for which ^ = 0), there occurs a 
merger of a stable and an unstable branch of the (0, 1) type solution. Thus, in other words, 
the (0, 1) solution is stable for values of K (and p) in between the points A and B, while 
it becomes unstable at B and the p-K graph corresponding to the unstable (m = 0, / = 1) 
solution continues up to C, the latter being located on the graph of p^{K) against K. As K 
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is made to increase beyond the point C, on the other hand, one obtains yet another solution, 
namely one of type (m = 1, / = 2). The latter is an unstable solution that continues up to 
K ^ oo, the variation of p with K being shown schematically in the figure. 

In summary, then, fig. [T] depicts the graphs of p^{K) and p^{K) as functions of K, with 
reference to which the solutions of various types characterized by the integers m, I can be 
described. Each synchronized solution originates on either of the two curves (OiCi, O1C2 
in the figure) and the variation of p with K for the solution corresponds to a curve such 
AP, ABC or CQ. While the curve AP corresponds to the (m = 0, / = 0) solution, other 
solutions of type (m,0), for the various other possible values of m (refer to formula f llObp ) 
also correspond to similar curves. The figure shows schematically another such curve A'P' 
where A' is the point at which a stable (m, 1) type solution (for which the p-K curve extends 
along A'B') transforms into one of type (m,0). As already mentioned, all these (m, 0) type 
solutions are of the nature of stable phase waves for large values of K, the phases of the 
oscillators on the ring for such a solution being given by eq. ( !T6|) . Meanwhile, solutions for 
various other values of / are also described in terms of their own p-K curves. 

The following section includes a detailed consideration of the solutions for a particular case 
with three oscillators on a ring where we describe all the possible synchronized solutions. 
In the process we identify the critical value {K^) of K (corresponding to the point (C) of 
backbending of the curve ACB in fig. [1]) for the onset of synchronization. For K slightly 
below Kc, the system exhibits intermittent chaos. 

VII. A RING WITH THREE OSCILLATORS 

We consider a ring of = 3 oscillators with frequencies ui,U2,uj3, for which we assume, for 
the sake of concreteness, Ui = 1, U2 = —3, and a;3 = 2, in arbitrary units, corresponding to 
which one has a; = 0, si = (by definition), S2 = —l{= y) , S3 = 2(= x), . 

Synchronized solutions for this ring can be denoted symbolically by indicating whether the 
successive phase differences ipii^i = 1,2,3) belong to the interval [— |, |] (which we indicate by 
the symbol 'r', corresponding to the right half circle) or to ^] (symbol T, corresponding 
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to the left half circle) and, in adition, by the value of m (the value of / is determined by 
the number of '1' symbols in the string characterizing the intervals to which the successive 
phase difference belong). For instance, a 'rrr(O)' solution corresponds to all the three phase 
differences being of type 'r' (hence implying / = 0), with the sum of the three being 0. 
Similarly, a 'rll(l)' solution corresponds to -01 being of type 'r' and ip2, i^s being of type 
'1' (hence 1 = 2), with the sum of the three being 2n. The two intervals mentioned above 
include common end-points since these correspond to the points of transition of solutions of 
one type into another, such as the points A and C in fig. [H 

Making use of this notation one can make up a list of the possible synchronized solutions 
for the system under consideration. The possible values of / being / = 0, 1,2,3, the corre- 
sponding values of m in accordance with the formula fllOap are as folios 



/ = : m = 0, 
/ = 1 : m = 0,1, 
/ = 2 : m = 1, 

/ = 3:m = l,2. (22) 

The above list then includes a total of twelve possible solutions which we label as rrr(O), 
rrl(O), rrl(l), rh(0), rh(l), Irr(O), Irr(l), rll(l), Irl(l), nr(l), 111(1), and 111(2). 

Considering, to begin with, the rrr(O) synchronized solution for the sake of concreteness, and 
following notations as in sec. IIIII one finds that, for any given K, there exists a p{= sin ■0i) 
corresponding to such a solution, provided the following equation is satisfied. 



sin ^{p) + sin ^{p + uy) + sin ^{p + ux) = 0, (23) 

Here we use the notation u = for the sake of brevity and note that the function sin~^ 
denotes the principal value of the inverse sine, which corresponds to an angle lying in the 
range [— f , f] (this was denoted by the symbol a in eq. ( pa|) . while the symbol /3 was used 
to denote the inverse sine lying in the range [f , the two function branches agree at the 
argument value -|-1 ). Recall that, once p is determined from the above equation for any given 
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value of K (i.e., equivalently, of u), the phase angle differences ijji = (pi — 4>i-i{i = 1,2,3) 
are obtained as, respectively, the three terms on its left hand side, corresponding to which 
one obtains a one-parameter family of solutions for the phase angles (pi themselves (see 
sec. II Vp . and a one-parameter family of synchronized solutions is thereby obtained from 
equations ([2]), (jl]). In the following, while speaking of a synchronized solution (such as the 
rrr(O) solution under consideration) we will actually mean the corresponding one-parameter 
family. 

As mentioned in sec Wl\ sjTichronized solutions do not exist for K less than a certain 
threshold value Kq, corresponding to the point of intersection of the graphs of p+{K) and 
]?_ (fC) (refer to eq. ( |T8l) which, with our choice of parameters, is given by Kq = 3,p+ = 
p_ = po(say) = — |. As K is made to increase beyond Kq, the rrr(O) solution first makes its 
appearance at 



2{2x^ + - 2xy) y'^ - xy - x^/2x^ - 2xy 

A = Ai say = =,P = Pi say = ^ , 24 

^ 2x-y+^2x^-2xy ^ ' 2x^ + y^-2xy ^ ' 

where the point {Ki,pi) lies on the graph of pj^{K) (note the slight change in notation 

compared to sec. IVII where the symbol K^^'> was used in place of Ki). This means that, 

for K = Ki, the equation ( l23l) is satisfied with p = pi, for which "03 = f, the latter being 

one end point of the range [— f , f ]• Since this is also an end point of the range [f , ^], it is 

evident that {Ki,pi) also satisfies the equation 



sin ^{p) + sin ^ {p + uy) + {ir — sin '^{p + ux)) = 0, (25) 

which is precisely the equation determining a rrl(O) solution. However, the p-K curve for 
the rrl(O) solution extends to the left of A (fig. [1]), i.e., in other words, as K is made to 
increase from Kq, the rrl(O) solution gets transformed into rrr(O) as K crosses Ki. This p-K 
curve, describing the history of the rrl(O) solution has a point of backbending (the point B 
in fig. [H, see also fig. [2] obtained numerically) where there occurs a change in stability of this 
solution so that the segment CB corresponds to an unstable rrl(O) solution. 

The eigenvalues of the stability matrix M are seen to be 

K I 

Ai = 0, A2,3 = Y ( - (?i + 92 + gs) ± Y (g? + ql + ql- qiq2 - ^293 - qsqi) ) , (26) 
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and thus the value (Kc) of K for the point B corresponds to the condition 

qiq2 + 92^3 + ^391 = 0. (27) 

For our choice of the parameters cui, U2, OJ^, the condition (1271) corresponds to 

= 3.04224 (approx). (28) 

Looking at fig. [H one arrives at the important cnclusion that this is the critical value of 
K below which there exists no synchronized solution for the system under consideration 
while, as K is made to increase past Kc, there appears a pair of rrl(O) solutions by tangent 
bifurcation, of which one is a stable synchronized solution while the other is an unstable 
one. For K slightly below this critical value {K^) the system exhibits intermittent chaos (see 
section IVml). 

As K is made to increase beyond K^, the stable branch of the p-K curve describing the 
rrl(O) solution terminates on the p^{K)-K curve (the point A on the branch OiCi in fig. [T]), 
and the rrl(O) solution transforms to a rrr(O) one described by its own p-K curve (branch 
AP in fig. [H As already mentioned, in the large-X regime this p-K curve tends to the limit 
= "^2 = "^3 = (refer to eq. (I15al) ) which thus corresponds to a phase-locked solution for 
the system under consideration. 

Meanwhile, for K beyond K^ the unstable branch of the rrl(O) solution terminates on the 
P- (K) curve (point C in fig. [1]) at 

( \ 2(3:^ + 2y^ - 2xy) -z^ + xy- y^ly^ - 2xy 

K = /sTs say = ,p = P2 say = " , 29 

x-2y + A/2y2 - 2xy + 2y^ - 2xy 

where the rrl(O) solution gets transformed into a rll(l) one, an unstable synchronized solution 
for the system. The latter continues, for larger values of up to X —t- 00, when it tends to 
an unstable synchronized solution with ^1 = 0, ■02 = '*/'3 = Ti'- 
Figures |21 E] show the numerically computed p-K curves for the rrr(O), rrl(O), rll(l)type so- 
lutions, and the llr(l), 111(1), Irl(l) and 111(2) type solutions respectively, for the system with 
parameter values mentioned above. One observes that the 111(2) solution gets transformed 
into a Irl(l) solution and, similarly, a 111(1) solution gets transformed into a llr(l) solution in 
a manner analogous to the way a rrr(O) solution changes over to a rrl(O) one. Among these, 
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the 111(1) and 111(2) solutions are in the nature of synchronized phase waves on the ring for 
large with V'l = "02 = ''/'s = ^ and '0i = "02 = "03 = ^ respectively. Moreover, there 
occur points of backbending for the p-K curves describing the llr(l) and Irl(l) solutions, 
analogous to a similar backbending for the rrl(O) solution as seen above. These points once 
again correspond to the condition qiq2 + q2(l3, + qzli = (refer to eq. (J271l ). but now one 
of the two eigenvalues A2,3 remains positive while the other passes through zero, as a result 
of which the solution under consideration remains unstable on either side of the point of 
backbending. 
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-0.2 
-0.25 

A 

-0.3 
-0.35 
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-0.45 

3 3.05 3.1 3.15 3.2 3.25 3.3 3.35 3.4 3.45 3.5 

— > A" 

FIG. 2: p-K curves for the rrr(O), rrl(O) and rll(l) synchronized solutions; the point of backbending 
for the rrl(O) curve corresponds to the point C in fig. [U where a stable and an unstable branch 
originate by a tangent bifurcation. 

While figures |2l [3] depict the history of seven of the twelve possible synchronized solutions 
mentioned above, the remaining five (Irr(O), Irr(l), rlr(O), rlr(l), and rrl(l) ) happen to be 
ruled out for the particular situation under consideration in this section, since in each case, 
the existence of a solution to the equation determining p as a function of K (a particular 
instance of the general equation' f l9ap ) may be seen to be inconsistent with the requirement 
-l<p<+lioiKQ<K< oo. 
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FIG. 3: p-K curves for the llr(l), 111(1), Irl(l), and 111(2) synchronized solutions; the point of 
backbending on either of the llr(l) and Irl(l) solutions is analogous to that on the rrl(O) curve of 
fig. [21 where one eigenvalue of the stability matrix changes sign. 

In other words, for a ring with = 3 oscillators with the parameter values (wi, Ci;2, cus) 
chosen as above, there exist only seven synchronized solutions, of which one exists within a 
finite range of the parameter K while the other six persist up to large values of K. Among 
the latter six, only one is a stable synchronized solution with (pi = (f)2 = 03 (phase-locking) 
while the others are all unstable. 



In a recent work on a closely related model 22| , Mehta and Kastner have worked out fixeed 



points of a certain functional (the lattice Landau gauge fixing functional) in the context of 
lattice gauge theory that correspond to the synchronized solutions of the Kuramoto ring con- 
sidered here, and their results also relate to the stability characteristics of the synchronized 



solutions. The extension to two dimensional systems has also been considered [23 1. 
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VIII. INTERMITTENT CHAOS 



Fig. m depicts the alteration of laminar and turbulent phases in the intermittennt chaos 



24j |: the time evolution of the phase angles in the laminar phase may be approximated by 



a discrete mapping) for K slightly below where we look at a trajectory initiated close 
to the point on the 3-torus with t/^i, 7/'2, V'a corresponding to the rrl(O) solution at the point 
of the tangent bifurcation (at K = K^), and where we make the choice 0i = from among 
the one-parameter family of solutions corresponding to these values of -02, '03 (this gives 
cj)^ = 0, 02 ~ -1.39002335, 03 ^ 0.33238947). 
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FIG. 4: Intermittent chaos for K slightly below the critical value Kc{^ 3.04224); alteration of 
laminar and turbulent phases as seen in a numerical computation of the evolution of the phase 
angle 6'i with K = K^- 0.00390625; the initial condition (6'i(t = 0) = (jjf^ {i = 1,2,3)) has been 
chosen as explained in the text; the laminar phase corresponds to values of 9i close to <j)f \ 0^*''* + ^, 
0(°) + f . 



As seen in fig. HI the laminar phase alternates between values of 6i close to (f)\ \ + , 
^(0) _|_ ^{0) ]^gjj^g ^YiQ initial value of 9i chosen close to (the deviation from was chosen 
with a random number generator). This is explained by the fact that 6*1 + ^2 + ^3 is a constant 
of motion for the system under consideration and , at the same time, a constant shift $ in 
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the three phase angles results in the same orbit, taking these two facts togethr, one finds 
that the constant shift $ can be ^ or ^ for any given initial condition. 

It may be noted that, near the tangent bifurcation, the system behaves effectively as a one 
dimensional one since, among the three angles 9i, 62, 63, one is trivially eliminated by the 
condition + 62 + O3 = ujt) (which, moreover, is zero for our choice of the parameters) 
following from equations f|Ta|) . fllb|) . and then a centre manifold reduction leads to a single 
effective angular variable (say, G), with the remaining, third variable being slaved to it (recall 
that, close to the tangent bifurcation, one of the three eigenvalues of the stability matrix is 
trivially zero while, among the remaining two, the relevant eigenvalue crosses through zero 
and the third eigenvalue remains negative). The equation of motion for the relevant angular 
variable G in the laminar phase can then be assumed to be of the standard form 

G = e + G^ {e = K,-K). (30) 

As is evident from a dimensional argument, the duration of the laminar phase is then 
expected to scale with e as e~^. As seen in fig. 0, the scaling exponent obtained numerically 
for the system under consideration is indeed close to —1. 

The above phenomenon of intermittent chaos is a particular instance of the general phe- 
nomenon of chaos in systems of coupled oscillators before the onset of synchronization (see. 



e.g. 



[l5| cited above). In sec. IIX CI below we include a few numerical results on our 
system where the intermittent chaos is entrained to a synchronized solution by means of an 
external periodic driving. 

IX. PERIODIC DRIVING: ENTRAINMENT OF SYNCHRONIZED OSCILLA- 
TIONS 

A. Entrainment: a perturbative approach 

We now include a periodic driving term of amplitude A and frequency Q in eq. (Ilaj) so as 
to obtain the system 
K 

Oi = Ui + —{sin{e^+i - Oi) + sm(^i_i - 6^)) + Asin{VLt -9^ + 6), (i = 1, 2, . . . , A^; A > 0) 

(31) 
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FIG. 5: Scaling of the duration (r) of the laminar phase; r scales with e = — K as e'^, where 
u (ps —0.97) is seen to be close to the theoretical value —1. 

where 6 denotes a constant determining the difference between the driving phase {Qt+6) and 
the free-running phase angle (i.e., the phase angle in the absence of driving) characterizing 
the ith oscillator (z = 1,2,...,A^) and where the ring boundary conditions are implied. 
Driven populations of oscillators with a global coupling between them have been considered 
in which shows a rich pattern of behaviour for such systems, involving a number of 
bifurcation routes. In the present paper, we look for entrained solutions where the phases 
of the oscillators on the ring all vary periodically with the driving period Q as 

= nt + (Pi + -fi. (32) 

Here 0j denotes the constant part of the phase of the ith oscillator (z = 1,2, ... ,A^) for 
a synchronized solution of the form ([2]) which we take as a reference, calling the latter a 
free-running solution. The driving frequency Q characterizing the entrained solution differs 
from the free-running frequency u hj A = Q — u, which we refer to as the detuning. 
Finally, 7^ {i = 1,2,..., A^) are the deviations of the constant parts of the phases from the 
corresponding free- running values (pi. 
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We assume that the driving amphtude and the detuning are sufficiently small, as a result of 
which there corresponds to the free-running solution an entrained solution that differs by a 
small extent from the former, and outline a perturbative procedure for obtaining the phase 
deviations 7^ for the latter. Substituting the proposed solution fl52]) in f lTa|) one obtains 

2 

sm(0i+i -(pi + 7i+i - 7i) + sin{(pi_i - 0, + - = — (fi - Ui - Asin{5 -(pi- %)). 

(33) 

This system of equation is consistent only if 

N 

^(f] -u^- Asin{5 -(pi- 7i)) = 0. (34) 
1=1 

Making use of this condition, an approximation scheme for working out the solution for 
sufficiently small A and S can be set up as follows. 

We ffist choose appropiate variables in place of 7^ (z = 1, . . . , A^) by referring to the lin- 
earization of the left hand side of (!33|) . On making use of eq. ([3]), this can be expressed in 
the form Mi 7, where Mi is precisely times the stability matrix M given by ffTTl) (with 
(pi {i = 1,...,N) corresponding to the synchronized solution under consideration), and 
where 7 stands for the column 7 = {71, 72, ■ ■ ■ , 7Ar}"'"- Since M possesses a zero eigenvalue, a 
solution of eq. f l33p by straightforward linearization and matrix inversion is not possible. If 
T be the orthogonal matrix diagonalizing M, then the normal co-ordinates for the linearized 
problem (say, pi, p2, • • • , Pn) are given by 

p = {pi,...,p;vr = r% (35) 

If the eigenvalues of M are denoted by (Ai = 0, A2, . . . , Xn) (i.e., those of Mi are -^Aj (z = 
1,2,..., A^)), then the ffist column of T, corresponding to the eigenvector for the eigenvalue 
zero, is made up of -^'s as all the entries, and one has 



1 ^ 

7X7 E 7. (36) 



. 1 

1=1 



where is a normalization factor. 
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The linearization of the left hand side of eq. (133|) would then give (refer to equa- 
tions (I6a|). fl6b|) ) 

Ap = T^C, (37) 

where C = {Ci, C2, • • • , Ca^}^, with 

Q = ^{A-Asin{6- (Pi (38) 
and A is the diagonal form of Mi, with -^Aj as the ith diagonal element. 

Noting that the first row of is just "^{l? ^, ■ ■ ■ A}^^ oiie finds that the first of the 
equations in (!37j) is an identity (in view of eq. (IMl) ). which is why one has to retain the 
nonlinear terms in C. Indeed, the choice of the new variables pi{i = 1,2, ... , N) is a useful 
one because, as we see below, it allows us to make a clear separation between large and 
small components in the column p. It turns out that pi is large compared to the remaining 
components and includes a uniform phase translation for all the oscillators in the ring. In 
contrast to the other components of p, pi is appropriately solved for from the consistency 
relation eq. flMl) . The latter can be written in the form (recall that Xlj*^* ~ Nu) 



N 



J]] (A - Asm(A - 0i - {Tp)i)) = (A = - w). (39) 



i=l 



The basic idea underlying our approximation scheme is that, for sufficiently small A and 
A, P2, ■ ■ ■ , Pn remain small, being obtained from the — 1 equations of ( 1371) with the first 
equation omitted, while pi, which is expected to be large compared to p2i • • • , Ptv, is to be 
solved for from eq. fl39|) . Accordingly, we first define 

7'°' = (rp)„.,^^.„.„, (40) 

i.e., expressing 7 in terms of the pj's, we retain only pi, ignoring the rest of the normal 
co-ordinates. Noting, however, that the first column of T is just -^{^^ ■ ■ ■ , 1}^, we obtain 
(for z = 1,2, ...,A^) 

7f ^ = (41) 



A^ 
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The required entrained solution is then obtained in the following steps. 



Substituting x = • • • ' Pu ^'i- pSi) . we solve for pi from eq. (l39|l . i.e. from 

N 

(A - - 0. - -^)) = 0. (42) 

Evidently, the solution for pi arrived at represents only the leading order approximation, 
which we denote as pf'\ In this leading order, p2, ■ ■ ■ , Pn are all zero. What we are interested 
in is the next order of approximation to pi as also to p2, ■ ■ ■ , Pn- 

One obtains p2, ■ ■ ■ , Pn from the second,. . .,A^th equations in eq. fl37|) where, in the right 
hand side, we substitute (^^^ for (, being obtained from eq. ( 138!) by replacing 7^ with 
7} = (i = 1, . . . , A^). In other words, one has 



where 



P. = f Ari(rV°0^, {^ = 2,3,...,N) (43a) 



Cr = -(A-Aszn(5-0.-;^)). (43b) 



Finally, denoting the correction to pi over p^^^ by r;, 

Pi=pS°^+r/, (44) 

one determines t] by referring to the consistency condition (l34ll . now written with the cor- 
rection terms included in = {Tp)i, i.e, with p written as {pf^ + ri,p2,. . . jPn}"^, where 
P2, . . ■ , Pn are determined from eq. (143 al) as above. Having obtained now a consistent ap- 
proximation for all components of p, one can finally work out 7 = Tp, and then the entrained 
solution fl32D. 



It may be noted that the effect of the periodic driving on a synchronized solution is two-fold, 
involving, first, the phase shifts = 1,2, .. . , A^) for the oscillators modifying the time- 
independent phases (pi and, secondly, the change in the frequency from the free-running 
frequency u to the driving frequency Q. Of these the phase shifts 7^ are made up of a 
dominant term 7^°^ which is the same for all the oscillators, along with correction terms 
determined by rj and p2, ■ ■ ■ , Pn as given by eq. (I43ap . However, a uniform phase shift for 
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all the oscillators does not imply any change in the orbit of the system of oscillators in 
its phase space (a A^-torus) since the orbit already represents a one-parameter family of 
solutions corresponding to various values of a uniform phase shift $. Thus, the relevant 
part of the effect of the driving is contained in the correction terms mentioned above, apart 
from the shift in the frequency. 



B. Frequency window for a given driving amplitude 

An immediate consequence of the consistency condition fl34l) is that the detuning A and the 
driving amplitude A have to satisfy 

|A| < A, (45) 

for an entrained solution (eq. f l32|) ) to exist. In other words, for a given value of the driving 
amplitude A, the detuning A has to be within the window —A < A < A. This condition 
requiring the detuning to be within a certain window depending on the driving amplitude 



is a common feature in experiments and models on entrainment (see, e.g., 25|). 



Fig. E] depicts an entrained solution obtained by the above perturbative approach by starting 
from a free-running stable rrr(O) solution with A^ = 3. Referring to sec. IIX El below, one finds 
that this is a stable entrained solution for the system under consideration. One observes 
from the figure that the perturbation approach works excellently when compared with the 
time evolution obtained by a numerical integration of the system under consideration. 



C. Entrainment of intermittent chaos 



Fig. [71 on the other hand, depicts the phenomenon of entrainment of intermittent chaos 
where a periodic driving results in an entrained solution for the system (with A^ = 3) for 
parameter values corresponding to intermittent chaos {K < K^). The driving frequency Q 
for the existance of such an entrained solution has to satisfy once again the condition f H5]) 
where now u{= Q — A) is to be interpreted as the mean of the oscillator frequencies rather 
than the free-running frequency of an entrained solution. 
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FIG. 6: Depicting a stable entrained solution obtained by perturbation from a free-running rrr(O) 
one; K = 10.8066978, 5 = 0,^ = 0.03125, A = 0.015625; the curve frawn with a solid line cor- 
responds to the solution obtained by the perturbation approach explained in sec. HXl while the 
crosses are points generated by numerical integration; the agreement between the two is evident. 

By contrast, fig. [8] depicts a situation where the amplitude and frequency of the driving term 
violates the condition fH5l) . with K slightly less than K^., in which one finds that entrainment 
does not take place, resulting instead in an irregular variation of the phase angles. 

D. Entrained solutions for large A 

While we have worked out the entrained solutions for small values of the driving amplitude 
(and, hence, of the detuning), there evidently exists entrained solutions for larger values of 
A as well. However, the solutions for large A do not necessarily have any correspondence 
with the free-running synchronized solutions obtained earlier. Thus, for A» K, one has, 
clS cl first approximation, the entrained solution 



where the time- independent phases 0i, which now have no reference to the free-running 




(46a) 
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FIG. 7: Entrainment of intermittent chaos; K = — 6, A, and A as in fig. [6l where A is now 
to be interpreted as the difference between 17 and the mean of the oscillator frequencies, rather than 
a free-running synchronized frequency; an initial condition is chosen close to 9i = (f>i mentioned in 
the first paragraph of sec. IVIII^ by numerically solving for 7^ {i = 1, 2, 3) in eq. ([33|) . 



phases, are given by 



0. = 5-[.zn-i](^L_^). (46b) 



Here the function [sin ^] differs from the principal value of the inverse sine in that it may 
correspond to either of the two branches a and /3 introduced in sec. IIIH i.e., to an angle in 



either of the two ranges [0, |) and 



1-2 ' 2 ^ 



The condition for the existence of such an entrained solution is now (see eq. (15b|) ) 

\n - uji\ = \Ai\ < A {i = 1,2,..., N), (47) 
which, at the same time, implies condition f H5|) . 



Evidently, the above solution is stable if all the 0j's are determined from eq. fl46b|) by 
choosing the branch a for [sm~^]. 

For small values of ^, entrained solutions can be obtained perturbatively from the above 
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FIG. 8: Irregular time variation of one of the three phase angles in the case of periodic driving on 
intermittent chaos {K slightly less than Kc] we choose K = — g^g) where the detuning is not 
small enough to satisfy the condition given by ()45p (as in fig. [TJ is to be interpreted as the mean 
of the oscillator frequencies rather than the frequency of a free-running synchronized solution) ; we 
choose A = A + for the purpose of illustration. 

solutions determined for A ^ oo. Numerically, one obtains a host of such entrained solutions 
even away from the perturbative regime ^ << 1. However, the problem of continuation to 
the regime ^ << 1 remains. 

E. Stability of the entrained solutions 

Having obtained an entrained solution for either ^ << lor-^ << 1, or even away from these 
perturbative regimes by numerical computation, one can address the question of stability of 
the solution, which turns out to be a trivial one since the linearized equations are of the 
form 




(48) 
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Here {i = 1,2, . . . , N) denote the deviations from the phases (pi of the entrained solution 
(eq. ( I46ap ) under consideration, and the matrix M' is given (for = 3, for instance,) by 

(-{qi + q2 + ^cos (f)[) " „ \ 



M' 



K 



(49) 



q2 -{q2 + ga + ^cos 02) gs 

?i ^3 -{q-i + qi + ^cos 03) J 

In this expression, stands for cos^cpi — 0i_i) and (f)[ for (0j — 5) (« = 1,2, .. . , A^). In 
the perturbative regime ^ << 1, one has to substitute for 0j the perturbed phase angles 
{(pi + 7j) where, in the latter expression, (pi stands for the phase angles of the free-running 
solution and 7^ for the perturbations considered in sec. HXl 

The fact that the stability matrix M' is a constant rather than a periodic one reduces the 
question of stability to a simple determination of the eigenvalues of M'. 




FIG. 9: Depicting a stable entrained solution obtained numerically for K = 5, A = 10, 5 = 0, = 
5; these values correspond to the system being away from both the perturbative regimes ^ << 1 



and ^ << 1. 



Fig. [9] depicts an entrained solution obtained numerically ioi K = 5, A = 10, Q = 5, S = 0, 
for which (H7|) is satisfied. The initial condition for the soilution was obtained by numerically 
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solving eq. (133|) for 0i + 7j {i = 1,2,3), starting from values given by eq. (??) (where the 
principal value of sin~^ was taken for all the three angles), and is a stable solution by virtue 
of all the eigenvalues of M' (eq. (H^ ) being negative. It may be noted that the values of the 
parameters chosen correspond to the system being away from both the perturbative regimes 
mentioned above. 
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